/* Copyright 2007-2008 dnAnalytics Project.
 *
 * Contributors to this file:
 * Patrick van der Velde
 * Marcus Cuda - http://marcus.cuda.net
 * Jurgen Van Gael 
 *
 * This file is part of dnAnalytics.  dnAnalytics is licensed under the 
 * Microsoft Public License. See License.txt for a complete copy of the
 * license.
 */

using System;
using System.Collections.Generic;
using System.Text;
using dnAnalytics.LinearAlgebra.Decomposition;
using dnAnalytics.Math;
using dnAnalytics.Properties;

namespace dnAnalytics.LinearAlgebra
{
    /// <summary>
    /// A matrix class.
    /// </summary>
    public abstract class Matrix : IFormattable, ICloneable, IEquatable<Matrix>
    {
        private int mColumns;
        private int mRows;

        /// <summary>
        /// Sets the number of rows and columns for the matrix.
        /// </summary>
        /// <param name="rows">The number of rows.</param>
        /// <param name="columns">The number of columns.</param>
        protected Matrix(int rows, int columns)
        {
            if (rows < 1)
            {
                throw new ArgumentException(Resources.NotPositive, "rows");
            }
            if (columns < 1)
            {
                throw new ArgumentException(Resources.NotPositive, "columns");
            }
            mRows = rows;
            mColumns = columns;
        }

        /// <summary>
        /// Set or get an element of the matrix.
        /// </summary>
        /// <param name="row">The row index.</param>
        /// <param name="column">The column index.</param>
        /// <returns>A specific element in the matrix.</returns>
        /// <exception cref="ArgumentOutOfRangeException">If either the <paramref name="row"/> or <paramref name="column"/> index is outside the 
        /// the dimensions of this matrix.</exception>
        public virtual double this[int row, int column]
        {
            get
            {
                RangeCheck(row, column);
                return ValueAt(row, column);
            }
            set
            {
                RangeCheck(row, column);
                ValueAt(row, column, value);
            }
        }

        /// <summary>
        /// The number of columns.
        /// </summary>
        public int Columns
        {
            get { return mColumns; }
            protected set { mColumns = value; }
        }

        /// <summary>
        /// The numbers of rows.
        /// </summary>        
        public int Rows
        {
            get { return mRows; }
            protected set { mRows = value; }
        }

        #region System.Object overrides

        /// <summary>
        /// Compares two matrices for equality by checking that corresponding matrix entries are exactly equal.
        /// </summary>
        public override bool Equals(object obj)
        {
            return Equals(obj as Matrix);
        }

        /// <summary>
        /// Compares two matrices for equality by checking that corresponding matrix entries are exactly equal.
        /// </summary>
        public bool Equals(Matrix other)
        {
            // Reject equality when the argument is null or has a different shape.
            if (other == null)
            {
                return false;
            }
            if (Columns != other.Columns || Rows != other.Rows)
            {
                return false;
            }

            // Accept if the argument is the same object as this.
            if (ReferenceEquals(this, other))
            {
                return true;
            }

            // If all else fails, perform elementwise comparison.
            for (int i = 0; i < Rows; i++)
            {
                for (int j = 0; j < Columns; j++)
                {
                    if (ValueAt(i, j) != other.ValueAt(i, j))
                    {
                        return false;
                    }
                }
            }

            return true;
        }

        /// <summary>
        /// Computes a hash for the matrix based on 25 entries using the standard .NET hash computation for double arrays.
        /// </summary>
        public override int GetHashCode()
        {
            int hashNum = System.Math.Max(Rows * Columns, 25);
            double[] hashBase = new double[hashNum];
            for (int i = 0; i < hashNum; i++)
            {
                int col = i % Columns;
                int row = (i - col) / Rows;
                hashBase[i] = this[row, col];
            }
            return hashBase.GetHashCode();
        }

        #endregion

        #region ICloneable Members

        object ICloneable.Clone()
        {
            return Clone();
        }

        #endregion

        #region IFormattable Members

        ///<summary>A string representation of this matrix.</summary>
        ///<param name="format">A format specification.</param>
        ///<param name="formatProvider">An <see cref="IFormatProvider"/> that supplies culture-specific formatting information.</param>
        ///<returns>The string representation of this matrix as comma seperated values.</returns>
        public virtual string ToString(string format, IFormatProvider formatProvider)
        {
            StringBuilder sb = new StringBuilder();
            for (int i = 0; i < Rows; i++)
            {
                for (int j = 0; j < Columns; j++)
                {
                    sb.Append(ValueAt(i, j).ToString(format, formatProvider));
                    if (j != Columns - 1)
                    {
                        sb.Append(", ");
                    }
                }
                if (i != Rows - 1)
                {
                    sb.Append(Environment.NewLine);
                }
            }
            return sb.ToString();
        }

        #endregion

        /// <summary>Calculates the L1 norm.</summary>
        /// <returns>The L1 norm of the matrix.</returns>
        public virtual double L1Norm()
        {
            double norm = 0;
            for (int j = 0; j < mColumns; j++)
            {
                double s = 0;
                for (int i = 0; i < mRows; i++)
                {
                    s += System.Math.Abs(ValueAt(i, j));
                }
                norm = System.Math.Max(norm, s);
            }
            return norm;
        }

        /// <summary>Calculates the L2 norm.</summary>
        /// <returns>The L2 norm of the matrix.</returns>   
        /// <remarks>For sparse matrices, the L2 norm is computed using a dense implementation of singular value decomposition. 
        /// In a later release, it will be replaced with a sparse implementation.</remarks>
        public virtual double L2Norm()
        {
            Svd svd = new Svd(this, false);
            return svd.Norm2();
        }

        /// <summary>Calculates the Frobenius norm of this matrix.</summary>
        /// <returns>The Frobenius norm of this matrix.</returns>
        public virtual double FrobeniusNorm()
        {
            Matrix aat = this*Transpose();
            double norm = 0;
            for (int i = 0; i < mRows; i++)
            {
                norm += System.Math.Abs(aat.ValueAt(i, i));
            }
            norm = System.Math.Sqrt(norm);

            return norm;
        }

        /// <summary>Calculates the infinity norm of this matrix.</summary>
        /// <returns>The infinity norm of this matrix.</returns>   
        public virtual double InfinityNorm()
        {
            double norm = 0;
            for (int i = 0; i < mRows; i++)
            {
                double s = 0;
                for (int j = 0; j < mColumns; j++)
                {
                    s += System.Math.Abs(ValueAt(i, j));
                }
                norm = System.Math.Max(norm, s);
            }
            return norm;
        }

        /// <summary>Calculates the condition number of this matrix.</summary>
        /// <returns>The condition number of the matrix.</returns>
        /// <remarks>The condition number is calculated using singular value decomposition.</remarks>
        /// <exception cref="ConvergenceFailedException">If the singular value decomposition
        /// cannot converge to a solution.</exception>
        public virtual double ConditionNumber()
        {
            Svd svd = new Svd(this, false);
            return svd.ConditionNumber();
        }

        /// <summary>
        /// Sets all values to zero.
        /// </summary>
        public virtual void Clear()
        {
            for (int i = 0; i < Rows; i++)
            {
                for (int j = 0; j < Columns; j++)
                {
                    ValueAt(i, j, 0);
                }
            }
        }

        /// <summary>
        /// Copies a column into a new <see cref="Vector"/>.
        /// </summary>
        /// <param name="index">The column to copy.</param>
        /// <returns>A <see cref="Vector"/> containing the copied elements.</returns>
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="index"/> is negative,
        /// or greater than or equal to the number of columns.</exception>
        public virtual Vector GetColumn(int index)
        {
            Vector result = CreateVector(Rows);
            GetColumn(index, 0, Rows, result);
            return result;
        }

        /// <summary>
        /// Copies a column into to the given <see cref="Vector"/>.
        /// </summary>
        /// <param name="index">The column to copy.</param>
        /// <param name="result">The <see cref="Vector"/> to copy the column into.</param>
        /// <exception cref="ArgumentNullException">If the result <see cref="Vector"/> is <see langword="null" />.</exception>
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="index"/> is negative,
        /// or greater than or equal to the number of columns.</exception>
        /// <exception cref="NotConformableException">If <b>this.Rows != result.Count</b>.</exception>
        public virtual void GetColumn(int index, Vector result)
        {
            GetColumn(index, 0, Rows, result);
        }

        /// <summary>
        /// Copies the requested column elements into a new <see cref="Vector"/>.
        /// </summary>
        /// <param name="columnIndex">The column to copy elements from.</param>
        /// <param name="rowIndex">The row to start copying from.</param>
        /// <param name="length">The number of elements to copy.</param>
        /// <returns>A <see cref="Vector"/> containing the requested elements.</returns>
        /// <exception cref="ArgumentOutOfRangeException">If:
        /// <list><item><paramref name="columnIndex"/> is negative,
        /// or greater than or equal to the number of columns.</item>
        /// <item><paramref name="rowIndex"/> is negative,
        /// or greater than or equal to the number of rows.</item>
        /// <item><c>(rowIndex + length) &gt;= Rows.</c></item></list>
        /// </exception>        
        /// <exception cref="ArgumentException">If <paramref name="length"/> is not positive.</exception>
        public virtual Vector GetColumn(int columnIndex, int rowIndex, int length)
        {
            Vector result = CreateVector(length);
            GetColumn(columnIndex, rowIndex, length, result);
            return result;
        }

        /// <summary>
        /// Copies the requested column elements into the given vector.
        /// </summary>
        /// <param name="columnIndex">The column to copy elements from.</param>
        /// <param name="rowIndex">The row to start copying from.</param>
        /// <param name="length">The number of elements to copy.</param>
        /// <param name="result">The <see cref="Vector"/> to copy the column into.</param>
        /// <exception cref="ArgumentNullException">If the result <see cref="Vector"/> is <see langword="null" />.</exception>
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="columnIndex"/> is negative,
        /// or greater than or equal to the number of columns.</exception>        
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="rowIndex"/> is negative,
        /// or greater than or equal to the number of rows.</exception>        
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="rowIndex"/> + <paramref name="length"/>  
        /// is greater than or equal to the number of rows.</exception>
        /// <exception cref="ArgumentException">If <paramref name="length"/> is not positive.</exception>
        /// <exception cref="NotConformableException">If <strong>result.Count &lt; length</strong>.</exception>
        public virtual void GetColumn(int columnIndex, int rowIndex, int length, Vector result)
        {
            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (columnIndex >= mColumns || columnIndex < 0)
            {
                throw new ArgumentOutOfRangeException("columnIndex");
            }
            if (rowIndex >= mRows || rowIndex < 0)
            {
                throw new ArgumentOutOfRangeException("rowIndex");
            }

            if (rowIndex + length > mRows)
            {
                throw new ArgumentOutOfRangeException("length");
            }

            if (length < 1)
            {
                throw new ArgumentException(Resources.NotPositive, "length");
            }

            if (result.Count < length)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }

            for (int i = rowIndex, j = 0; i < rowIndex + length; i++, j++)
            {
                result[j] = ValueAt(i, columnIndex);
            }
        }

        /// <summary>
        /// Creates a new matrix and inserts the given column at the given index.
        /// </summary>
        /// <param name="columnIndex">The index of where to insert the column.</param>
        /// <param name="column">The column to insert.</param>
        /// <returns>A new matrix with the inserted column.</returns>
        /// <exception cref="ArgumentNullException">If <paramref name="column "/> is <see langword="null" />. </exception>
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="columnIndex"/> is &lt; zero or &gt; the number of columns.</exception>
        /// <exception cref="NotConformableException">If the size of <paramref name="column"/> != the number of rows.</exception>
        public virtual Matrix InsertColumn(int columnIndex, Vector column)
        {
            if (column == null)
            {
                throw new ArgumentNullException("column");
            }
            if (columnIndex < 0 || columnIndex > mColumns)
            {
                throw new ArgumentOutOfRangeException("columnIndex");
            }
            if (column.Count != mRows)
            {
                throw new NotConformableException("row", Resources.ParameterNotConformable);
            }
            Matrix result = CreateMatrix(mRows, mColumns + 1);
            for( int index = 0; index< columnIndex; index++)
            {
                result.SetColumn(index, GetColumn(index));
            }
            result.SetColumn(columnIndex, column);

            for( int index = columnIndex + 1; index < mColumns+1; index++)
            {
                result.SetColumn(index, GetColumn(index-1));
            }
            return result;
        }

        /// <summary>
        /// Returns an <see cref="IEnumerator{T}"/> that enumerates the requested matrix columns.
        /// </summary>
        /// <param name="index">The column to start enumerating over.</param>
        /// <param name="length">The number of columns to enumerating over.</param>
        /// <returns>An <see cref="IEnumerator{T}"/> that enumerates over requested matrix columns.</returns>
        /// <seealso cref="IEnumerator{T}"/>
        /// <exception cref="ArgumentOutOfRangeException">If:
        /// <list><item><paramref name="index"/> is negative,
        /// or greater than or equal to the number of columns.</item>        
        /// <item><c>(index + length) &gt;= Columns.</c></item></list>
        /// </exception>   
        /// <exception cref="ArgumentException">If <paramref name="length"/> is not positive.</exception>     
        public virtual IEnumerable<KeyValuePair<int, Vector>> GetColumnEnumerator(int index, int length)
        {
            if (index >= mColumns || index < 0)
            {
                throw new ArgumentOutOfRangeException("index");
            }
            if (index + length > mColumns)
            {
                throw new ArgumentOutOfRangeException("length");
            }
            if (length < 1)
            {
                throw new ArgumentException(Resources.NotPositive, "length");
            }
            int maxi = index + length;
            for (int i = index; i < maxi; i++)
            {
                yield return new KeyValuePair<int, Vector>(i, GetColumn(i));
            }
        }

        /// <summary>
        /// Returns an <see cref="IEnumerator{T}"/> that enumerates over the matrix columns.
        /// </summary>
        /// <returns>An <see cref="IEnumerator{T}"/> that enumerates over the matrix columns</returns>
        /// <seealso cref="IEnumerator{T}"/>
        public virtual IEnumerable<KeyValuePair<int, Vector>> GetColumnEnumerator()
        {
            for (int i = 0; i < mColumns; i++)
            {
                yield return new KeyValuePair<int, Vector>(i, GetColumn(i));
            }
        }

        /// <summary>
        /// Copies a row into an <see cref="Vector"/>.
        /// </summary>
        /// <param name="index">The row to copy.</param>
        /// <returns>A <see cref="Vector"/> containing the copied elements.</returns>
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="index"/> is negative,
        /// or greater than or equal to the number of rows.</exception>
        public virtual Vector GetRow(int index)
        {
            Vector ret = CreateVector(mColumns);
            GetRow(index, 0, Columns, ret);
            return ret;
        }

        /// <summary>
        /// Copies a row into to the given <see cref="Vector"/>.
        /// </summary>
        /// <param name="index">The row to copy.</param>
        /// <param name="result">The <see cref="Vector"/> to copy the row into.</param>
        /// <exception cref="ArgumentNullException">If the result vector is <see langword="null" />.</exception>
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="index"/> is negative,
        /// or greater than or equal to the number of rows.</exception>
        /// <exception cref="NotConformableException">If <b>this.Columns != result.Count</b>.</exception>
        public virtual void GetRow(int index, Vector result)
        {
            GetRow(index, 0, Columns, result);
        }

        /// <summary>
        /// Copies the requested row elements into a new <see cref="Vector"/>.
        /// </summary>
        /// <param name="rowIndex">The row to copy elements from.</param>
        /// <param name="columnIndex">The column to start copying from.</param>
        /// <param name="length">The number of elements to copy.</param>
        /// <returns>A <see cref="Vector"/> containing the requested elements.</returns>
        /// <exception cref="ArgumentOutOfRangeException">If:
        /// <list><item><paramref name="rowIndex"/> is negative,
        /// or greater than or equal to the number of rows.</item>
        /// <item><paramref name="columnIndex"/> is negative,
        /// or greater than or equal to the number of columns.</item>
        /// <item><c>(columnIndex + length) &gt;= Columns.</c></item></list></exception>        
        /// <exception cref="ArgumentException">If <paramref name="length"/> is not positive.</exception>
        public virtual Vector GetRow(int rowIndex, int columnIndex, int length)
        {
            Vector ret = CreateVector(length);
            GetRow(rowIndex, columnIndex, length, ret);
            return ret;
        }

        /// <summary>
        /// Copies the requested row elements into a new <see cref="Vector"/>.
        /// </summary>
        /// <param name="rowIndex">The row to copy elements from.</param>
        /// <param name="columnIndex">The column to start copying from.</param>
        /// <param name="length">The number of elements to copy.</param>
        /// <param name="result">The <see cref="Vector"/> to copy the column into.</param>
        /// <exception cref="ArgumentNullException">If the result <see cref="Vector"/> is <see langword="null" />.</exception>
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="rowIndex"/> is negative,
        /// or greater than or equal to the number of columns.</exception>        
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="columnIndex"/> is negative,
        /// or greater than or equal to the number of rows.</exception>        
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="columnIndex"/> + <paramref name="length"/>  
        /// is greater than or equal to the number of rows.</exception>
        /// <exception cref="ArgumentException">If <paramref name="length"/> is not positive.</exception>
        /// <exception cref="NotConformableException">If <strong>result.Count &lt; length</strong>.</exception>
        public virtual void GetRow(int rowIndex, int columnIndex, int length, Vector result)
        {
            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (rowIndex >= mRows || rowIndex < 0)
            {
                throw new ArgumentOutOfRangeException("rowIndex");
            }
            if (columnIndex >= mColumns || columnIndex < 0)
            {
                throw new ArgumentOutOfRangeException("columnIndex");
            }

            if (columnIndex + length > mColumns)
            {
                throw new ArgumentOutOfRangeException("length");
            }

            if (length < 1)
            {
                throw new ArgumentException(Resources.NotPositive, "length");
            }

            if (result.Count < length)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }

            for (int i = columnIndex, j = 0; i < columnIndex + length; i++, j++)
            {
                result[j] = ValueAt(rowIndex, i);
            }
        }

        /// <summary>
        /// Creates a new matrix and inserts the given row at the given index.
        /// </summary>
        /// <param name="rowIndex">The index of where to insert the row.</param>
        /// <param name="row">The row to insert.</param>
        /// <returns>A new matrix with the inserted column.</returns>
        /// <exception cref="ArgumentNullException">If <paramref name="row"/> is <see langword="null" />. </exception>
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="rowIndex"/> is &lt; zero or &gt; the number of rows.</exception>
        /// <exception cref="NotConformableException">If the size of <paramref name="row"/> != the number of columns.</exception>
        public virtual Matrix InsertRow(int rowIndex, Vector row)
        {
            if( row == null)
            {
                throw new ArgumentNullException("row");
            }
            if( rowIndex < 0 || rowIndex > mRows)
            {
                throw new ArgumentOutOfRangeException("rowIndex");
            }
            if( row.Count != mColumns)
            {
                throw new NotConformableException("row", Resources.ParameterNotConformable);
            }

            Matrix result = CreateMatrix(mRows+1, mColumns);
            for (int index = 0; index < rowIndex; index++)
            {
                result.SetRow(index, GetRow(index));
            }
            result.SetRow(rowIndex, row);

            for (int index = rowIndex + 1; index < mRows + 1; index++)
            {
                result.SetRow(index, GetRow(index - 1));
            }
            return result;
        }


        /// <summary>
        /// Returns an <see cref="IEnumerator{T}"/> that enumerates the requested matrix rows.
        /// </summary>
        /// <param name="index">The row to start enumerating over.</param>
        /// <param name="length">The number of rows to enumerating over.</param>
        /// <returns>An <see cref="IEnumerator{T}"/> that enumerates over requested matrix rows.</returns>
        /// <seealso cref="IEnumerator{T}"/>
        /// <exception cref="ArgumentOutOfRangeException">If:
        /// <list><item><paramref name="index"/> is negative,
        /// or greater than or equal to the number of rows.</item>        
        /// <item><c>(index + length) &gt;= Rows.</c></item></list></exception>        
        /// <exception cref="ArgumentException">If <paramref name="length"/> is not positive.</exception>     
        public virtual IEnumerable<KeyValuePair<int, Vector>> GetRowEnumerator(int index, int length)
        {
            if (index >= mRows || index < 0)
            {
                throw new ArgumentOutOfRangeException("index");
            }
            if (index + length > mRows)
            {
                throw new ArgumentOutOfRangeException("length");
            }
            if (length < 1)
            {
                throw new ArgumentException(Resources.NotPositive, "length");
            }
            int maxi = index + length;
            for (int i = index; i < maxi; i++)
            {
                yield return new KeyValuePair<int, Vector>(i, GetRow(i));
            }
        }

        /// <summary>
        /// Returns an <see cref="IEnumerator{T}"/> that enumerates over the matrix rows.
        /// </summary>
        /// <returns>An <see cref="IEnumerator{T}"/> that enumerates over the matrix rows</returns>
        /// <seealso cref="IEnumerator{T}"/>        
        public virtual IEnumerable<KeyValuePair<int, Vector>> GetRowEnumerator()
        {
            for (int i = 0; i < mRows; i++)
            {
                yield return new KeyValuePair<int, Vector>(i, GetRow(i));
            }
        }

        /// <summary>
        /// Copies the elements of this matrix to the given matrix.
        /// </summary>
        /// <param name="target">The matrix to copy values into.</param>
        /// <exception cref="ArgumentNullException">If target is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If target is not the same size as this matrix.</exception>
        public virtual void CopyTo(Matrix target)
        {
            if (target == null)
            {
                throw new ArgumentNullException("target");
            }

            if (Rows != target.Rows || Columns != target.Columns)
            {
                throw new NotConformableException("target", Resources.ParameterNotConformable);
            }

            foreach (KeyValuePair<int, Vector> column in GetColumnEnumerator())
            {
                foreach (KeyValuePair<int, double> element in column.Value)
                {
                    target.ValueAt(element.Key, column.Key, element.Value);
                }
            }
        }

        /// <summary>Computes the determinant of this matrix.</summary>
        /// <returns>The determinant of this matrix.</returns>
        /// <exception cref="MatrixNotSquareException">If the matrix is not square.</exception>
        /// <remarks>The determinant is calculated using LU decomposition.</remarks>
        public virtual double Determinant()
        {
            if (Rows != Columns)
            {
                throw new MatrixNotSquareException();
            }
            LU lu = new LU(this);
            return lu.Determinant();
        }

        /// <summary>Computes the trace of this matrix.</summary>
        /// <returns>The trace of this matrix.</returns>
        /// <exception cref="MatrixNotSquareException">If the matrix is not square.</exception>
        public virtual double Trace()
        {
            if (Rows != Columns)
            {
                throw new MatrixNotSquareException();
            }

            double t = 0.0;
            for (int i = 0; i < Rows; i++)
            {
                t += this[i, i];
            }

            return t;
        }

        /// <summary>
        /// Returns the elements of the diagonal in a <see cref="Vector"/>.
        /// </summary>
        /// <returns>The elements of the diagonal.</returns>
        /// <remarks>For non-square matrices, the method returns Min(Rows, Columns) elements where
        /// i == j (i is the the row index, and j is the column index).</remarks>
        public virtual Vector GetDiagonal()
        {
            int min = System.Math.Min(mRows, mColumns);
            Vector ret = CreateVector(min);
            for (int i = 0; i < min; i++)
            {
                ret[i] = ValueAt(i, i);
            }
            return ret;
        }

        /// <summary>
        /// Returns the transpose of this matrix.
        /// </summary>        
        /// <returns>The transpose of this matrix.</returns>
        public virtual Matrix Transpose()
        {
            Matrix ret = CreateMatrix(mColumns, mRows);
            for (int j = 0; j < mColumns; j++)
            {
                for (int i = 0; i < mRows; i++)
                {
                    ret.ValueAt(j, i, ValueAt(i, j));
                }
            }
            return ret;
        }

        /// <summary>
        /// Returns the inversls
        /// rman e of this matrix. The inverse is calculated using LU decomposition.
        /// </summary>
        /// <returns>The inverse of this matrix.</returns>
        /// <exception cref="MatrixNotSquareException">If this matrix is not square.</exception>
        public virtual Matrix Inverse()
        {
            if (mRows != mColumns)
            {
                throw new MatrixNotSquareException();
            }
            return new LU(this).Inverse();
        }

        /// <summary>
        /// Copies the values of the given <see cref="Vector"/> to the specified column.
        /// </summary>
        /// <param name="index">The column to copy the values to.</param>
        /// <param name="source">The vector to copy the values from.</param>
        /// <exception cref="ArgumentNullException">If <paramref name="source"/> is <see langword="null" />.</exception>        
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="index"/> is less than zero,
        /// or greater than or equal to the number of columns.</exception>
        /// <exception cref="NotConformableException">If the size of <paramref name="source"/> does not
        /// equal the number of rows of this <strong>Matrix</strong>.</exception>
        public virtual void SetColumn(int index, Vector source)
        {
            if (index < 0 || index >= mColumns)
            {
                throw new ArgumentOutOfRangeException("index");
            }

            if (source == null)
            {
                throw new ArgumentNullException("source");
            }

            if (source.Count != mRows)
            {
                throw new NotConformableException("source", Resources.ParameterNotConformable);
            }

            for (int i = 0; i < mRows; i++)
            {
                ValueAt(i, index, source[i]);
            }
        }

        /// <summary>
        /// Copies the values of the given array to the specified column.
        /// </summary>
        /// <param name="index">The column to copy the values to.</param>
        /// <param name="source">The array to copy the values from.</param>
        /// <exception cref="ArgumentNullException">If <paramref name="source"/> is <see langword="null" />.</exception>        
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="index"/> is less than zero,
        /// or greater than or equal to the number of columns.</exception>
        /// <exception cref="NotConformableException">If the size of <paramref name="source"/> does not
        /// equal the number of rows of this <strong>Matrix</strong>.</exception>
        /// <exception cref="NotConformableException">If the size of <paramref name="source"/> does not
        /// equal the number of rows of this <strong>Matrix</strong>.</exception>
        public virtual void SetColumn(int index, double[] source)
        {
            if (index < 0 || index >= mColumns)
            {
                throw new ArgumentOutOfRangeException("index");
            }

            if (source == null)
            {
                throw new ArgumentNullException("source");
            }

            if (source.Length != mRows)
            {
                throw new NotConformableException("source", Resources.ParameterNotConformable);
            }

            for (int i = 0; i < mRows; i++)
            {
                ValueAt(i, index, source[i]);
            }
        }

        /// <summary>
        /// Copies the values of the given <see cref="Vector"/> to the diagonal.
        /// </summary>
        /// <param name="source">The vector to copy the values from. The length of the vector should be
        /// Min(Rows, Columns).</param>
        /// <exception cref="ArgumentNullException">If <paramref name="source"/> is <see langword="null" />.</exception>   
        /// <exception cref="NotConformableException">If the length of <paramref name="source"/> does not
        /// equal Min(Rows, Columns).</exception>
        /// <remarks>For non-square matrices, the elements of <paramref name="source"/> are copied to
        /// this[i,i].</remarks>
        public virtual void SetDiagonal(Vector source)
        {
            if (source == null)
            {
                throw new ArgumentNullException("source");
            }

            int min = System.Math.Min(Rows, Columns);
            if (source.Count != min)
            {
                throw new NotConformableException("source", Resources.ParameterNotConformable);
            }
            for (int i = 0; i < min; i++)
            {
                ValueAt(i, i, source[i]);
            }
        }

        /// <summary>
        /// Copies the values of the given array to the diagonal.
        /// </summary>
        /// <param name="source">The array to copy the values from. The length of the array should be
        /// Min(Rows, Columns).</param>
        /// <exception cref="ArgumentNullException">if <paramref name="source"/> is <see langword="null" />.</exception>        
        /// <exception cref="NotConformableException">If the length of <paramref name="source"/> does not
        /// equal Min(Rows, Columns).</exception>
        /// <remarks>For non-square matrices, the elements of <paramref name="source"/> are copied to
        /// this[i,i].</remarks>
        public virtual void SetDiagonal(double[] source)
        {
            if (source == null)
            {
                throw new ArgumentNullException("source");
            }

            int min = System.Math.Min(Rows, Columns);
            if (source.Length != min)
            {
                throw new NotConformableException("source", Resources.ParameterNotConformable);
            }
            for (int i = 0; i < min; i++)
            {
                ValueAt(i, i, source[i]);
            }
        }

        /// <summary>
        /// Copies the values of the given <see cref="Vector"/> to the specified row.
        /// </summary>
        /// <param name="index">The row to copy the values to.</param>
        /// <param name="source">The vector to copy the values from.</param>
        /// <exception cref="ArgumentNullException">If <paramref name="source"/> is <see langword="null" />.</exception>            
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="index"/> is less than zero,
        /// or greater than or equal to the number of rows.</exception>
        /// <exception cref="NotConformableException">If the size of <paramref name="source"/> does not
        /// equal the number of columns of this <strong>Matrix</strong>.</exception>
        public virtual void SetRow(int index, Vector source)
        {
            if (index < 0 || index >= mRows)
            {
                throw new ArgumentOutOfRangeException("index");
            }

            if (source == null)
            {
                throw new ArgumentNullException("source");
            }

            if (source.Count != mColumns)
            {
                throw new NotConformableException("source", Resources.ParameterNotConformable);
            }

            for (int j = 0; j < mColumns; j++)
            {
                ValueAt(index, j, source[j]);
            }
        }

        /// <summary>
        /// Copies the values of the given array to the specified row.
        /// </summary>
        /// <param name="index">The row to copy the values to.</param>
        /// <param name="source">The array to copy the values from.</param>
        /// <exception cref="ArgumentNullException">If <paramref name="source"/> is <see langword="null" />.</exception>  
        /// <exception cref="ArgumentOutOfRangeException">If <paramref name="index"/> is less than zero,
        /// or greater than or equal to the number of rows.</exception>
        /// <exception cref="NotConformableException">If the size of <paramref name="source"/> does not
        /// equal the number of columns of this <strong>Matrix</strong>.</exception>
        public virtual void SetRow(int index, double[] source)
        {
            if (index < 0 || index >= mRows)
            {
                throw new ArgumentOutOfRangeException("index");
            }

            if (source == null)
            {
                throw new ArgumentNullException("source");
            }

            if (source.Length != mColumns)
            {
                throw new NotConformableException("source", Resources.ParameterNotConformable);
            }

            for (int j = 0; j < mColumns; j++)
            {
                ValueAt(index, j, source[j]);
            }
        }

        /// <summary>
        /// Returns a new matrix containing the lower triangle of this matrix.
        /// </summary>
        /// <returns>The lower triangle of this matrix.</returns>        
        public virtual Matrix GetLowerTriangle()
        {
            Matrix ret = CreateMatrix(mRows, mColumns);
            for (int j = 0; j < mColumns; j++)
            {
                for (int i = 0; i < mRows; i++)
                {
                    if (i >= j)
                    {
                        ret.ValueAt(i, j, ValueAt(i, j));
                    }
                }
            }
            return ret;
        }

        /// <summary>
        /// Puts the the lower triangle of this matrix into the result matrix.
        /// </summary>
        /// <param name="result">Where to store the the lower triangle.</param>
        /// <exception cref="ArgumentNullException">If <paramref name="result"/> is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not the same as this matrix.</exception>
        public virtual void GetLowerTriangle(Matrix result)
        {
            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (result.Rows != Rows || result.Columns != Columns)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }
            for (int j = 0; j < mColumns; j++)
            {
                for (int i = 0; i < mRows; i++)
                {
                    if (i >= j)
                    {
                        result.ValueAt(i, j, ValueAt(i, j));
                    }
                    else
                    {
                        result.ValueAt(i, j, 0);
                    }
                }
            }
        }

        /// <summary>
        /// Returns a new matrix containing the lower triangle of this matrix. The new matrix
        /// does not contain the diagonal elements of this matrix.
        /// </summary>
        /// <returns>The lower triangle of this matrix.</returns>
        public virtual Matrix GetStrictlyLowerTriangle()
        {
            Matrix ret = CreateMatrix(mRows, mColumns);
            for (int j = 0; j < mColumns; j++)
            {
                for (int i = 0; i < mRows; i++)
                {
                    if (i > j)
                    {
                        ret.ValueAt(i, j, ValueAt(i, j));
                    }
                }
            }
            return ret;
        }

        /// <summary>
        /// Puts the the strictly lower triangle of this matrix into the result matrix.
        /// </summary>
        /// <param name="result">Where to store the the lower triangle.</param>
        /// <exception cref="ArgumentNullException">If <paramref name="result"/> is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not the same as this matrix.</exception>
        public virtual void GetStrictlyLowerTriangle(Matrix result)
        {
            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (result.Rows != Rows || result.Columns != Columns)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }
            for (int j = 0; j < mColumns; j++)
            {
                for (int i = 0; i < mRows; i++)
                {
                    if (i > j)
                    {
                        result.ValueAt(i, j, ValueAt(i, j));
                    }
                    else
                    {
                        result.ValueAt(i, j, 0);
                    }
                }
            }
        }

        /// <summary>
        /// Returns a new matrix containing the upper triangle of this matrix.
        /// </summary>
        /// <returns>The upper triangle of this matrix.</returns>   
        public virtual Matrix GetUpperTriangle()
        {
            Matrix ret = CreateMatrix(mRows, mColumns);
            for (int j = 0; j < mColumns; j++)
            {
                for (int i = 0; i < mRows; i++)
                {
                    if (i <= j)
                    {
                        ret.ValueAt(i, j, ValueAt(i, j));
                    }
                }
            }
            return ret;
        }

        /// <summary>
        /// Puts the the upper triangle of this matrix into the result matrix.
        /// </summary>
        /// <param name="result">Where to store the the lower triangle.</param>
        /// <exception cref="ArgumentNullException">If <paramref name="result"/> is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not the same as this matrix.</exception>
        public virtual void GetUpperTriangle(Matrix result)
        {
            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (result.Rows != Rows || result.Columns != Columns)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }
            for (int j = 0; j < mColumns; j++)
            {
                for (int i = 0; i < mRows; i++)
                {
                    if (i <= j)
                    {
                        result.ValueAt(i, j, ValueAt(i, j));
                    }
                    else
                    {
                        result.ValueAt(i, j, 0);
                    }
                }
            }
        }

        /// <summary>
        /// Returns a new matrix containing the upper triangle of this matrix. The new matrix
        /// does not contain the diagonal elements of this matrix.
        /// </summary>
        /// <returns>The upper triangle of this matrix.</returns>
        public virtual Matrix GetStrictlyUpperTriangle()
        {
            Matrix ret = CreateMatrix(mRows, mColumns);
            for (int j = 0; j < mColumns; j++)
            {
                for (int i = 0; i < mRows; i++)
                {
                    if (i < j)
                    {
                        ret.ValueAt(i, j, ValueAt(i, j));
                    }
                }
            }
            return ret;
        }

        /// <summary>
        /// Puts the the strictly upper triangle of this matrix into the result matrix.
        /// </summary>
        /// <param name="result">Where to store the the lower triangle.</param>
        /// <exception cref="ArgumentNullException">If <paramref name="result"/> is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not the same as this matrix.</exception>
        public virtual void GetStrictlyUpperTriangle(Matrix result)
        {
            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (result.Rows != Rows || result.Columns != Columns)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }
            for (int j = 0; j < mColumns; j++)
            {
                for (int i = 0; i < mRows; i++)
                {
                    if (i < j)
                    {
                        result.ValueAt(i, j, ValueAt(i, j));
                    }
                    else
                    {
                        result.ValueAt(i, j, 0);
                    }
                }
            }
        }

        /// <summary>
        /// Creates a matrix that contains the values from the requested sub-matrix.
        /// </summary>
        /// <param name="rowIndex">The row to start copying from.</param>
        /// <param name="rowLength">The number of rows to copy. Must be positive.</param>
        /// <param name="columnIndex">The column to start copying from.</param>
        /// <param name="columnLength">The number of columns to copy. Must be positive.</param>
        /// <exception cref="ArgumentOutOfRangeException">If: <list><item><paramref name="rowIndex"/> is
        /// is negative, or greater than or equal to the number of rows.</item>
        /// <item><paramref name="columnIndex"/> is is negative, or greater than or equal to the number 
        /// of columns.</item>
        /// <item><c>(columnIndex + columnLength) &gt;= Columns</c></item>
        /// <item><c>(rowIndex + rowLength) &gt;= Rows</c></item></list></exception>        
        /// <exception cref="ArgumentException">If <paramref name="rowLength"/> or <paramref name="columnLength"/>
        /// is not positive.</exception>
        public virtual Matrix GetSubMatrix(int rowIndex, int rowLength, int columnIndex, int columnLength)
        {
            if (rowIndex >= mRows || rowIndex < 0)
            {
                throw new ArgumentOutOfRangeException("rowIndex");
            }

            if (columnIndex >= mColumns || columnIndex < 0)
            {
                throw new ArgumentOutOfRangeException("columnIndex");
            }

            if (rowLength < 1)
            {
                throw new ArgumentException(Resources.NotPositive, "rowLength");
            }
            if (columnLength < 1)
            {
                throw new ArgumentException(Resources.NotPositive, "columnLength");
            }
            int colMax = columnIndex + columnLength;
            int rowMax = rowIndex + rowLength;

            if (rowMax > mRows)
            {
                throw new ArgumentOutOfRangeException("rowLength");
            }

            if (colMax > mColumns)
            {
                throw new ArgumentOutOfRangeException("columnLength");
            }

            Matrix result = CreateMatrix(rowLength, columnLength);
            for (int j = columnIndex, jj = 0; j < colMax; j++, jj++)
            {
                for (int i = rowIndex, ii = 0; i < rowMax; i++, ii++)
                {
                    result.ValueAt(ii, jj, ValueAt(i, j));
                }
            }
            return result;
        }

        /// <summary>
        /// Returns this matrix as a multidimensional array.
        /// </summary>
        /// <returns>A multidimensional containing the values of this matrix.</returns>        
        public virtual double[,] ToArray()
        {
            double[,] ret = new double[mRows,mColumns];
            for (int j = 0; j < mColumns; j++)
            {
                for (int i = 0; i < mRows; i++)
                {
                    ret[i, j] = ValueAt(i, j);
                }
            }
            return ret;
        }

        /// <summary>
        /// Returns the matrix's elements as an array with the data laid out column-wise.
        /// </summary>
        /// <example><pre>
        /// 1, 2, 3
        /// 4, 5, 6  will be returned as  1, 4, 7, 2, 5, 8, 3, 6, 9
        /// 7, 8, 9
        /// </pre></example>
        /// <returns>An array containing the matrix's elements.</returns>
        public virtual double[] ToColumnWiseArray()
        {
            double[] ret = new double[Rows*Columns];
            foreach (KeyValuePair<int, Vector> column in GetColumnEnumerator())
            {
                int jIndex = column.Key*Rows;
                foreach (KeyValuePair<int, double> element in column.Value)
                {
                    ret[jIndex + element.Key] = element.Value;
                }
            }

            return ret;
        }

        /// <summary>
        /// Returns the matrix's elements as an array with the data laid row-wise.
        /// </summary>
        /// <example><pre>
        /// 1, 2, 3
        /// 4, 5, 6  will be returned as  1, 2, 3, 4, 5, 6, 7, 8, 9
        /// 7, 8, 9
        /// </pre></example>
        /// <returns>An array containing the matrix's elements.</returns>
        public virtual double[] ToRowWiseArray()
        {
            double[] ret = new double[Rows*Columns];

            foreach (KeyValuePair<int, Vector> row in GetRowEnumerator())
            {
                int iIndex = row.Key*Columns;
                foreach (KeyValuePair<int, double> element in row.Value)
                {
                    ret[iIndex + element.Key] = element.Value;
                }
            }
            return ret;
        }

        /// <summary>
        /// Returns a clone of this matrix.
        /// </summary>
        /// <returns>A clone of this matrix.</returns>
        /// <remarks>Added as an alternative to the unary addition operator.</remarks>
        public virtual Matrix Plus()
        {
            return Clone();
        }

        /// <summary>
        /// Adds a scalar to each element in the matrix.
        /// </summary>
        /// <param name="scalar">The scalar to add.</param>        
        public virtual void Add(double scalar)
        {
            if (FloatingPointUtilities.EqualsWithinSignificantFigures(0.0, scalar, 15))
            {
                return;
            }

            for (int i = 0; i < Rows; i++)
            {
                for (int j = 0; j < Columns; j++)
                {
                    ValueAt(i, j, ValueAt(i, j) + scalar);
                }
            }
        }

        /// <summary>
        /// Adds a scalar to each element in the matrix and places the results into the result matrix.
        /// </summary>
        /// <param name="scalar">The scalar to add.</param>
        /// <param name="result">The results of the addition.</param>
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not the same as this matrix.</exception>
        public virtual void Add(double scalar, Matrix result)
        {
            if (result == null)
            {
                throw new ArgumentNullException("result");
            }
            if (result.Rows != Rows || result.Columns != Columns)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }
            CopyTo(result);
            result.Add(scalar);
        }

        /// <summary>
        /// Adds another matrix to this matrix.
        /// </summary>
        /// <param name="other">The matrix to add to this matrix.</param>
        /// <exception cref="ArgumentNullException">If the other matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If the two matrices don't have the same dimensions.</exception>
        public virtual void Add(Matrix other)
        {
            if (other == null)
            {
                throw new ArgumentNullException("other");
            }

            if (other.Rows != Rows || other.Columns != Columns)
            {
                throw new NotConformableException("other", Resources.ParameterNotConformable);
            }

            for (int i = 0; i < Rows; i++)
            {
                for (int j = 0; j < Columns; j++)
                {
                    double tmp = ValueAt(i, j);
                    tmp += other.ValueAt(i, j);
                    ValueAt(i, j, tmp);
                }
            }
        }

        /// <summary>
        /// Adds another matrix to this matrix and places the results into the result matrix.
        /// </summary>
        /// <param name="other">The matrix to add to.</param>
        /// <param name="result">The results of the addition.</param>
        /// <exception cref="ArgumentNullException">If the other matrix is <see langword="null" />.</exception>
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If the two matrices don't have the same dimensions.</exception>
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not the same as this matrix.</exception>
        public virtual void Add(Matrix other, Matrix result)
        {
            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (result.Rows != Rows || result.Columns != Columns)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }
            CopyTo(result);
            result.Add(other);
        }

        /// <summary>
        /// Negates each element of this matrix.
        /// </summary>        
        public virtual void Negate()
        {
            Multiply(-1);
        }

        /// <summary>
        /// Negate each element of this matrix and place the results into to result matrix.
        /// </summary>
        /// <param name="result">The result of the negation.</param>
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">if the result matrix's dimensions are not the same as this matrix.</exception>
        public virtual void Negate(Matrix result)
        {
            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (result.Rows != Rows || result.Columns != Columns)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }
            CopyTo(result);
            result.Negate();
        }

        /// <summary>
        /// Subtracts a scalar from each value in the matrix.
        /// </summary>
        /// <param name="scalar">The scalar to subtract.</param>
        public virtual void Subtract(double scalar)
        {
            Add(-scalar);
        }

        /// <summary>
        /// Subtracts a scalar from each in the matrix and place the results into the result matrix.
        /// </summary>
        /// <param name="scalar">The scalar to subtract.</param>
        /// <param name="result">The result of the subtraction.</param>
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not the same as this matrix.</exception>
        public virtual void Subtract(double scalar, Matrix result)
        {
            if (result == null)
            {
                throw new ArgumentNullException("result");
            }
            if (result.Rows != Rows || result.Columns != Columns)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }
            CopyTo(result);
            result.Subtract(scalar);
        }

        /// <summary>
        /// Subtracts another matrix from this matrix.
        /// </summary>
        /// <param name="other">The matrix to subtract.</param>
        /// <exception cref="ArgumentNullException">If the other matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If the two matrices don't have the same dimensions.</exception>
        public virtual void Subtract(Matrix other)
        {
            if (other == null)
            {
                throw new ArgumentNullException("other");
            }

            if (other.Rows != Rows || other.Columns != Columns)
            {
                throw new NotConformableException("other", Resources.ParameterNotConformable);
            }

            for (int i = 0; i < Rows; i++)
            {
                for (int j = 0; j < Columns; j++)
                {
                    double tmp = ValueAt(i, j);
                    tmp -= other.ValueAt(i, j);
                    ValueAt(i, j, tmp);
                }
            }
        }

        /// <summary>
        /// Subtracts another matrix for this matrix and places the results into the result matrix.
        /// </summary>
        /// <param name="other">The matrix to subtract.</param>
        /// <param name="result">The results of the subtraction.</param>
        /// <exception cref="ArgumentNullException">If the other matrix is <see langword="null" />.</exception>
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If the two matrices don't have the same dimensions.</exception>
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not the same as this matrix.</exception>
        public virtual void Subtract(Matrix other, Matrix result)
        {
            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (result.Rows != Rows || result.Columns != Columns)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }
            CopyTo(result);
            result.Subtract(other);
        }

        /// <summary>
        /// Multiplies each element of this matrix with a scalar.
        /// </summary>
        /// <param name="scalar">The scalar to multiply with.</param>
        public virtual void Multiply(double scalar)
        {
            if (FloatingPointUtilities.EqualsWithinSignificantFigures(1.0, scalar, 15))
            {
                return;
            }
            for (int i = 0; i < Rows; i++)
            {
                for (int j = 0; j < Columns; j++)
                {
                    ValueAt(i, j, ValueAt(i, j)*scalar);
                }
            }
        }

        /// <summary>
        /// Pointwise multiplies this matrix with another matrix.
        /// </summary>
        /// <param name="other">The matrix to pointwise multiply with this one.</param>
        /// <exception cref="ArgumentNullException">If the other matrix is <see langword="null" />.</exception> 
        /// <exception cref="NotConformableException">If this matrix and <paramref name="other"/> are not the same size.</exception>
        /// <returns>A new matrix that is the pointwise multiplication of this matrix and <paramref name="other"/>.</returns>
        public virtual Matrix PointwiseMultiply(Matrix other)
        {
            if (other == null)
            {
                throw new ArgumentNullException("other");
            }
            if (Columns != other.Columns || Rows != other.Rows)
            {
                throw new NotConformableException("other", Resources.ParameterNotConformable);
            }

            Matrix result = CreateMatrix(Rows, Columns);
            PointwiseMultiply(other, result);
            return result;
        }

        /// <summary>
        /// Pointwise multiplies this matrix with another matrix and stores the result into the result matrix.
        /// </summary>
        /// <param name="other">The matrix to pointwise multiply with this one.</param>
        /// <param name="result">The matrix to store the result of the pointwise multiplication.</param>
        /// <exception cref="ArgumentNullException">If the other matrix is <see langword="null" />.</exception> 
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception> 
        /// <exception cref="NotConformableException">If this matrix and <paramref name="other"/> are not the same size.</exception>
        /// <exception cref="NotConformableException">If this matrix and <paramref name="result"/> are not the same size.</exception>
        public virtual void PointwiseMultiply(Matrix other, Matrix result)
        {
            if (other == null)
            {
                throw new ArgumentNullException("other");
            }
            
            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (Columns != other.Columns || Rows != other.Rows)
            {
                throw new NotConformableException("result");
            }

            if (Columns != result.Columns || Rows != result.Rows)
            {
                throw new NotConformableException("result");
            }

            for (int i = 0; i < Rows; i++)
            {
                for (int j = 0; j < Columns; j++)
                {
                    result.ValueAt(i, j, ValueAt(i, j) * other.ValueAt(i, j));
                }
            }
        }

        /// <summary>
        /// Multiplies each element of the matrix by a scalar and places results into the result matrix.
        /// </summary>
        /// <param name="scalar">The scalar to multiply the matrix with.</param>
        /// <param name="result">The matrix to multiply.</param>
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception> 
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not the same as this matrix.</exception>
        public virtual void Multiply(double scalar, Matrix result)
        {
            if (result == null)
            {
                throw new ArgumentNullException("result");
            }
            if (result.Rows != Rows || result.Columns != Columns)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }

            CopyTo(result);
            result.Multiply(scalar);
        }

        /// <summary>
        /// Multiplies this matrix with another matrix and returns the result.
        /// </summary>
        /// <param name="other">The matrix to multiply with.</param>
        /// <exception cref="NotConformableException">If <strong>this.Columns != other.Rows</strong>.</exception>        
        /// <exception cref="ArgumentNullException">If the other matrix is <see langword="null" />.</exception>
        public virtual Matrix Multiply(Matrix other)
        {
            if (other == null)
            {
                throw new ArgumentNullException("other");
            }

            if (Columns != other.Rows)
            {
                throw new NotConformableException("other", Resources.ParameterNotConformable);
            }
            Matrix result = CreateMatrix(Rows, other.Columns);
            Multiply(other, result);
            return result;
        }

        /// <summary>
        /// Multiplies this matrix with another matrix and places the results into the result matrix.
        /// </summary>
        /// <param name="other">The matrix to multiply with.</param>
        /// <param name="result">The result of the multiplication.</param>
        /// <exception cref="ArgumentNullException">If the other matrix is <see langword="null" />.</exception>
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If <strong>this.Columns != other.Rows</strong>.</exception>
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not the this.Rows x other.Columns.</exception>
        public virtual void Multiply(Matrix other, Matrix result)
        {
            if (other == null)
            {
                throw new ArgumentNullException("other");
            }

            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (Columns != other.Rows)
            {
                throw new NotConformableException("other", Resources.ParameterNotConformable);
            }

            if (result.Rows != Rows || result.Columns != other.Columns)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }

            for (int j = 0; j != other.Columns; j++)
            {
                for (int i = 0; i != Rows; i++)
                {
                    double s = 0;
                    for (int l = 0; l != Columns; l++)
                    {
                        s += ValueAt(i, l)*other.ValueAt(l, j);
                    }
                    result.ValueAt(i, j, s);
                }
            }
        }

        /// <summary>
        /// Multiplies this matrix by a vector and returns the result.
        /// </summary>
        /// <param name="rightSide">The vector to multiply with.</param>
        /// <returns>The result of the multiplication.</returns>
        /// <exception cref="ArgumentNullException">If rightSide is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If <c>this.Columns != rightSide.Size</c>.</exception>
        public virtual Vector Multiply(Vector rightSide)
        {
            Vector ret = CreateVector(Rows);
            Multiply(rightSide, ret);
            return ret;
        }

        /// <summary>
        /// Multiplies this matrix with a vector and places the results into the result matrix.
        /// </summary>
        /// <param name="rightSide">The vector to multiply with.</param>
        /// <param name="result">The result of the multiplication.</param>
        /// <exception cref="ArgumentNullException">If rightSide is <see langword="null" />.</exception>
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If <strong>result.Count != this.Rows</strong>.</exception>
        /// <exception cref="NotConformableException">If <strong>this.Columns != rightSide.Count</strong>.</exception>
        public virtual void Multiply(Vector rightSide, Vector result)
        {
            if (rightSide == null)
            {
                throw new ArgumentNullException("rightSide");
            }

            if (Columns != rightSide.Count)
            {
                throw new NotConformableException("rightSide", Resources.ParameterNotConformable);
            }

            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (Rows != result.Count)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }

            for (int i = 0; i != Rows; i++)
            {
                double s = 0;
                for (int j = 0; j != Columns; j++)
                {
                    s += ValueAt(i, j)*rightSide[j];
                }
                result[i] = s;
            }
        }

        /// <summary>
        /// Multiplies this matrix by a vector and returns the result. The operation is vector * matrix.
        /// </summary>
        /// <param name="leftSide">The vector to multiply with.</param>
        /// <returns>The result of the multiplication.</returns>
        /// <exception cref="ArgumentNullException">If leftSide is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If <strong>this.Rows != leftSide.Size</strong>.</exception>
        public virtual Vector LeftMultiply(Vector leftSide)
        {
            Vector ret = CreateVector(Columns);
            LeftMultiply(leftSide, ret);
            return ret;
        }

        /// <summary>
        /// Multiplies this matrix with a vector and places the results into the result matrix. The operation is vector * matrix.
        /// </summary>
        /// <param name="leftSide">The vector to multiply with.</param>
        /// <param name="result">The result of the multiplication.</param>
        /// <exception cref="ArgumentNullException">If leftSide is <see langword="null" />.</exception>
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If <strong>result.Count != this.Columns</strong>.</exception>
        /// <exception cref="NotConformableException">If <strong>this.Rows != leftSide.Count</strong>.</exception>
        public virtual void LeftMultiply(Vector leftSide, Vector result)
        {
            if (leftSide == null)
            {
                throw new ArgumentNullException("leftSide");
            }

            if (Rows != leftSide.Count)
            {
                throw new NotConformableException("leftSide", Resources.ParameterNotConformable);
            }

            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (Columns != result.Count)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }

            for (int j = 0; j != Columns; j++)
            {
                double s = 0;
                for (int i = 0; i != leftSide.Count; i++)
                {
                    s += leftSide[i]*ValueAt(i, j);
                }
                result[j] = s;
            }
        }

        /// <summary>
        /// Divides each element of the matrix by a scalar.
        /// </summary>
        /// <param name="scalar">The scalar to divide by.</param>
        public virtual void Divide(double scalar)
        {
            Multiply(1/scalar);
        }

        /// <summary>
        /// Divides this matrix by a scalar and places the results into the results matrix.
        /// </summary>
        /// <param name="scalar">The scalar to divide the matrix by.</param>
        /// <param name="result">The result of the division.</param>
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not the same as this matrix.</exception>
        public virtual void Divide(double scalar, Matrix result)
        {
            Multiply(1/scalar, result);
        }

        /// <summary>
        /// Stacks this matrix on top of the given matrix and returns the combined matrix.
        /// </summary>
        /// <param name="lower">The matrix to stack this matrix upon.</param>
        /// <returns>The combined matrix.</returns>
        /// <exception cref="ArgumentNullException">If lower is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If <strong>upper.Columns != lower.Columns</strong>.</exception>
        public virtual Matrix Stack(Matrix lower)
        {
            if (lower == null)
            {
                throw new ArgumentNullException("lower");
            }

            if (lower.Columns != Columns)
            {
                throw new NotConformableException("lower", Resources.ParameterNotConformable);
            }

            Matrix result = CreateMatrix(Rows + lower.Rows, Columns);
            Stack(lower, result);
            return result;
        }

        /// <summary>
        /// Stacks this matrix on top of the given matrix and places the result into the result matrix.
        /// </summary>
        /// <param name="lower">The matrix to stack this matrix upon.</param>
        /// <param name="result">The combined matrix.</param>
        /// <exception cref="ArgumentNullException">If lower is <see langword="null" />.</exception>
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If <strong>upper.Columns != lower.Columns</strong>.</exception>
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not (this.Rows + lower.rows) x this.Columns.</exception>
        public virtual void Stack(Matrix lower, Matrix result)
        {
            if (lower == null)
            {
                throw new ArgumentNullException("lower");
            }
            if (lower.Columns != Columns)
            {
                throw new NotConformableException("lower", Resources.ParameterNotConformable);
            }

            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (result.Rows != (Rows + lower.Rows) || result.Columns != Columns)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }


            for (int j = 0; j < Columns; j++)
            {
                for (int i = 0; i < Rows; i++)
                {
                    result.ValueAt(i, j, ValueAt(i, j));
                }
            }

            for (int j = 0; j < lower.Columns; j++)
            {
                for (int i = 0; i < lower.Rows; i++)
                {
                    result.ValueAt(i + Rows, j, lower.ValueAt(i, j));
                }
            }
        }

        /// <summary>
        /// Concatenates this matrix with the given matrix.
        /// </summary>
        /// <param name="right">The matrix to concatenate.</param>
        /// <returns>The combined matrix.</returns>
        /// <exception cref="ArgumentNullException">If right is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If <strong>this.Rows != right.Rows.</strong></exception>
        public virtual Matrix Append(Matrix right)
        {
            if (right == null)
            {
                throw new ArgumentNullException("right");
            }

            if (right.Rows != Rows)
            {
                throw new NotConformableException("right", Resources.ParameterNotConformable);
            }

            Matrix result = CreateMatrix(Rows, Columns + right.Columns);
            Append(right, result);
            return result;
        }

        /// <summary>
        /// Concatenates this matrix with the given matrix and places the result into the result matrix.
        /// </summary>
        /// <param name="right">The matrix to concatenate.</param>
        /// <param name="result">The combined matrix.</param>
        /// <exception cref="ArgumentNullException">If right is <see langword="null" />.</exception>
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If <strong>this.Rows != right.Rows</strong>.</exception>
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not this.Rows x (this.Columns + right.Columns).</exception>
        public virtual void Append(Matrix right, Matrix result)
        {
            if (right == null)
            {
                throw new ArgumentNullException("right");
            }

            if (right.Rows != Rows)
            {
                throw new NotConformableException("right", Resources.ParameterNotConformable);
            }

            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (result.Columns != (Columns + right.Columns) || result.Rows != Rows)
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }

            for (int j = 0; j < Columns; j++)
            {
                for (int i = 0; i < Rows; i++)
                {
                    result.ValueAt(i, j, ValueAt(i, j));
                }
            }

            for (int j = 0; j < right.Columns; j++)
            {
                for (int i = 0; i < right.Rows; i++)
                {
                    result.ValueAt(i, j + Columns, right.ValueAt(i, j));
                }
            }
        }

        /// <summary>
        /// Diagonally stacks his matrix on top of the given matrix. The new matrix is a M-by-N matrix, 
        /// where M = this.Rows + lower.Rows and N = this.Columns + lower.Columns.
        /// The values of off the off diagonal matrices/blocks are set to zero.
        /// </summary>
        /// <param name="lower">The lower, right matrix.</param>
        /// <exception cref="ArgumentNullException">If lower is <see langword="null" />.</exception>
        /// <returns>the combined matrix</returns>
        public virtual Matrix DiagonalStack(Matrix lower)
        {
            if (lower == null)
            {
                throw new ArgumentNullException("lower");
            }

            Matrix result = CreateMatrix(Rows + lower.Rows, Columns + lower.Columns);
            DiagonalStack(lower, result);
            return result;
        }

        /// <summary>
        /// Diagonally stacks his matrix on top of the given matrix and places the combined matrix into the result matrix.
        /// </summary>
        /// <param name="lower">The lower, right matrix.</param>
        /// <param name="result">The combined matrix</param>
        /// <exception cref="ArgumentNullException">If lower is <see langword="null" />.</exception>
        /// <exception cref="ArgumentNullException">If the result matrix is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If the result matrix's dimensions are not (this.Rows + lower.rows) x (this.Columns + lower.Columns).</exception>
        public virtual void DiagonalStack(Matrix lower, Matrix result)
        {
            if (lower == null)
            {
                throw new ArgumentNullException("lower");
            }

            if (result == null)
            {
                throw new ArgumentNullException("result");
            }

            if (result.Rows != (Rows + lower.Rows) || result.Columns != (Columns + lower.Columns))
            {
                throw new NotConformableException("result", Resources.ParameterNotConformable);
            }

            for (int j = 0; j < Columns; j++)
            {
                for (int i = 0; i < Rows; i++)
                {
                    result.ValueAt(i, j, ValueAt(i, j));
                }
            }

            for (int j = 0; j < lower.Columns; j++)
            {
                for (int i = 0; i < lower.Rows; i++)
                {
                    result.ValueAt(i + Rows, j + Columns, lower.ValueAt(i, j));
                }
            }
        }

        ///<summary>A string representation of this matrix.</summary>
        ///<returns>The string representation of this matrix as comma seperated values.</returns>
        public override string ToString()
        {
            return ToString(null, null);
        }

        ///<summary>A string representation of this matrix.</summary>
        ///<param name="format">A format specification.</param>
        ///<returns>The string representation of this matrix as comma seperated values.</returns>
        public virtual string ToString(string format)
        {
            return ToString(format, null);
        }

        ///<summary>A string representation of this matrix.</summary>
        ///<param name="formatProvider">An <see cref="IFormatProvider"/> that supplies culture-specific formatting information.</param>
        ///<returns>The string representation of this matrix as comma seperated values.</returns>
        public virtual string ToString(IFormatProvider formatProvider)
        {
            return ToString(null, formatProvider);
        }

        /// <summary>
        /// Returns a deep-copy clone of the matrix.
        /// </summary>
        /// <returns>A deep-copy clone of the matrix.</returns>
        public virtual Matrix Clone()
        {
            Matrix result = CreateMatrix(Rows, Columns);
            CopyTo(result);
            return result;
        }

        /// <summary>
        /// Checks to see if the indices are valid, if not it throws an ArgumentOutOfRangeException.
        /// </summary>
        private void RangeCheck(int row, int column)
        {
            if (row < 0 || row >= mRows)
            {
                throw new ArgumentOutOfRangeException("row");
            }
            if (column < 0 || column >= mColumns)
            {
                throw new ArgumentOutOfRangeException("column");
            }
        }

        /// <summary>
        /// Retrieves the value at the given indices. <strong>This method is not ranged checked.</strong> If the
        /// indices are outside the dimensions of the matrix, it may return an incorrect value or throw
        /// an <see cref="IndexOutOfRangeException"/>.
        /// </summary>
        /// <param name="row">The row index.</param>
        /// <param name="column">The column index.</param>
        /// <returns>The value at the given indices.</returns>
        /// <exception cref="IndexOutOfRangeException">Depending on the implementation, an <see cref="IndexOutOfRangeException"/>
        /// may be thrown if one of the indices is outside the dimensions of the matrix.</exception>
        protected internal abstract double ValueAt(int row, int column);

        /// <summary>
        /// Set the value at the given indices to the given value. <strong>This method is not ranged checked.</strong> If the
        /// indices are outside the dimensions of the matrix, it may set an incorrect element or throw
        /// an <see cref="IndexOutOfRangeException"/>.
        /// </summary>
        /// <param name="row">The row index.</param>
        /// <param name="column">The column index.</param>
        /// <param name="value">The value to set.</param>
        /// <exception cref="IndexOutOfRangeException">Depending on the implementation, an <see cref="IndexOutOfRangeException"/>
        /// may be thrown if one of the indices is outside the dimensions of the matrix.</exception>
        protected internal abstract void ValueAt(int row, int column, double value);

        /// <summary>
        /// Creates a <strong>Matrix</strong> for the given number of rows and columns.
        /// </summary>
        /// <param name="numberOfRows">The number of rows.</param>
        /// <param name="numberOfColumns">The number of columns.</param>
        /// <returns>A <strong>Matrix</strong> with the given dimensions.</returns>
        protected internal abstract Matrix CreateMatrix(int numberOfRows, int numberOfColumns);

        /// <summary>
        /// Creates a <strong>Vector</strong> with a the given dimension.
        /// </summary>
        /// <param name="size">The size of the vector.</param>
        /// <returns>A <strong>Vector</strong> with the given dimension.</returns>
        protected internal abstract Vector CreateVector(int size);

        /// <summary>
        /// Returns a <strong>Matrix</strong> containing the same values of rightSide. 
        /// </summary>
        /// <remarks>This method is included for completeness.</remarks>
        /// <param name="rightSide">The matrix to get the values from.</param>
        /// <returns>A matrix containing a the same values as <paramref name="rightSide"/>.</returns>
        /// <exception cref="ArgumentNullException">If <paramref name="rightSide"/> is <see langword="null" />.</exception>
        public static Matrix operator +(Matrix rightSide)
        {
            if (rightSide == null)
            {
                throw new ArgumentNullException("rightSide");
            }
            return rightSide.Clone();
        }

        /// <summary>
        /// Adds two matrices together and returns the results.
        /// </summary>
        /// <param name="leftSide">One of the matrices to add.</param>
        /// <param name="rightSide">One of the matrices to add.</param>
        /// <returns>The result of the addition.</returns>
        /// <exception cref="NotConformableException">If <paramref name="leftSide"/> and <paramref name="rightSide"/> don't have the same dimensions.</exception>
        /// <exception cref="ArgumentNullException">If <paramref name="leftSide"/> or <paramref name="rightSide"/> is <see langword="null" />.</exception>
        public static Matrix operator +(Matrix leftSide, Matrix rightSide)
        {
            if (rightSide == null)
            {
                throw new ArgumentNullException("rightSide");
            }
            if (leftSide == null)
            {
                throw new ArgumentNullException("leftSide");
            }

            if (leftSide.Rows != rightSide.Rows || leftSide.Columns != rightSide.Columns)
            {
                throw new NotConformableException(Resources.ParametersNotConformable);
            }

            Matrix ret = leftSide.Clone();
            ret.Add(rightSide);
            return ret;
        }

        /// <summary>
        /// Negates each element of the matrix.
        /// </summary>
        /// <param name="rightSide">The matrix to negate.</param>
        /// <returns>A matrix containing the negated values.</returns>
        /// <exception cref="ArgumentNullException">If <paramref name="rightSide"/> is <see langword="null" />.</exception>
        public static Matrix operator -(Matrix rightSide)
        {
            if (rightSide == null)
            {
                throw new ArgumentNullException("rightSide");
            }
            Matrix ret = rightSide.Clone();
            ret.Negate();
            return ret;
        }

        /// <summary>
        /// Subtracts two matrices and returns the results.
        /// </summary>
        /// <param name="leftSide">The matrix to subtract from.</param>
        /// <param name="rightSide">The matrix to subtract.</param>
        /// <returns>The result of the subtraction.</returns>
        /// <exception cref="NotConformableException">If <paramref name="leftSide"/> and <paramref name="rightSide"/> don't have the same dimensions.</exception>
        /// <exception cref="ArgumentNullException">If <paramref name="leftSide"/> or <paramref name="rightSide"/> is <see langword="null" />.</exception>
        public static Matrix operator -(Matrix leftSide, Matrix rightSide)
        {
            if (rightSide == null)
            {
                throw new ArgumentNullException("rightSide");
            }
            if (leftSide == null)
            {
                throw new ArgumentNullException("leftSide");
            }

            if (leftSide.Rows != rightSide.Rows || leftSide.Columns != rightSide.Columns)
            {
                throw new NotConformableException(Resources.ParametersNotConformable);
            }

            Matrix ret = leftSide.Clone();
            ret.Subtract(rightSide);
            return ret;
        }

        /// <summary>
        /// Divides a <strong>Matrix</strong> by a constant and returns the result.
        /// </summary>
        /// <param name="leftSide">The matrix to divide.</param>
        /// <param name="rightSide">The constant to divide the matrix by.</param>
        /// <returns>The result of the division.</returns>
        /// <exception cref="ArgumentNullException">If <paramref name="leftSide"/> is <see langword="null" />.</exception>
        public static Matrix operator /(Matrix leftSide, double rightSide)
        {
            if (leftSide == null)
            {
                throw new ArgumentNullException("leftSide");
            }
            Matrix ret = leftSide.Clone();
            ret.Divide(rightSide);
            return ret;
        }

        /// <summary>
        /// Multiplies a <strong>Matrix</strong> by a constant and returns the result.
        /// </summary>
        /// <param name="leftSide">The matrix to multiply.</param>
        /// <param name="rightSide">The constant to multiply the matrix by.</param>
        /// <returns>The result of the multiplication.</returns>
        /// <exception cref="ArgumentNullException">If <paramref name="leftSide"/> is <see langword="null" />.</exception>
        public static Matrix operator *(Matrix leftSide, double rightSide)
        {
            if (leftSide == null)
            {
                throw new ArgumentNullException("leftSide");
            }
            Matrix ret = leftSide.Clone();
            ret.Multiply(rightSide);
            return ret;
        }

        /// <summary>
        /// Multiplies a <strong>Matrix</strong> by a constant and returns the result.
        /// </summary>
        /// <param name="leftSide">The matrix to multiply.</param>
        /// <param name="rightSide">The constant to multiply the matrix by.</param>
        /// <returns>The result of the multiplication.</returns>
        /// <exception cref="ArgumentNullException">If <paramref name="rightSide"/> is <see langword="null" />.</exception>
        public static Matrix operator *(double leftSide, Matrix rightSide)
        {
            if (rightSide == null)
            {
                throw new ArgumentNullException("rightSide");
            }
            Matrix ret = rightSide.Clone();
            ret.Multiply(leftSide);
            return ret;
        }

        /// <summary>
        /// Multiplies two matrices.
        /// </summary>
        /// <param name="leftSide">One of the matrices to multiply.</param>
        /// <param name="rightSide">One of the matrices to multiply.</param>
        /// <returns>The result of multiplication.</returns>
        /// <exception cref="ArgumentNullException">If <paramref name="leftSide"/> or <paramref name="rightSide"/> is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If the dimensions of <paramref name="leftSide"/> or <paramref name="rightSide"/> don't conform.</exception>
        public static Matrix operator *(Matrix leftSide, Matrix rightSide)
        {
            if (leftSide == null)
            {
                throw new ArgumentNullException("leftSide");
            }

            if (rightSide == null)
            {
                throw new ArgumentNullException("rightSide");
            }

            if (leftSide.Columns != rightSide.Rows)
            {
                throw new NotConformableException(Resources.ParametersNotConformable);
            }
            return leftSide.Multiply(rightSide);
        }

        /// <summary>
        /// Multiplies a <strong>Matrix</strong> and a <see cref="Vector"/>.
        /// </summary>
        /// <param name="leftSide">The matrix to multiply.</param>
        /// <param name="rightSide">The vector to multiply.</param>
        /// <returns>The result of multiplication.</returns>
        /// <exception cref="ArgumentNullException">If <paramref name="leftSide"/> or <paramref name="rightSide"/> is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If <strong>leftSide.Columns ! = rightSide.Count</strong>.</exception>
        public static Vector operator *(Matrix leftSide, Vector rightSide)
        {
            if (leftSide == null)
            {
                throw new ArgumentNullException("leftSide");
            }

            if (rightSide == null)
            {
                throw new ArgumentNullException("rightSide");
            }

            if (leftSide.Columns != rightSide.Count)
            {
                throw new NotConformableException(Resources.ParametersNotConformable);
            }

            return leftSide.Multiply(rightSide);
        }

        /// <summary>
        /// Multiplies a <see cref="Vector"/> and a <strong>Matrix</strong>.
        /// </summary>
        /// <param name="leftSide">The vector to multiply.</param>
        /// <param name="rightSide">The matrix to multiply.</param>
        /// <returns>The result of multiplication.</returns>
        /// <exception cref="ArgumentNullException">If <paramref name="leftSide"/> or <paramref name="rightSide"/> is <see langword="null" />.</exception>
        /// <exception cref="NotConformableException">If <strong>leftSide.Count != rightSide.Rows</strong>.</exception>
        public static Vector operator *(Vector leftSide, Matrix rightSide)
        {
            if (leftSide == null)
            {
                throw new ArgumentNullException("leftSide");
            }

            if (rightSide == null)
            {
                throw new ArgumentNullException("rightSide");
            }

            if (leftSide.Count != rightSide.Rows)
            {
                throw new NotConformableException(Resources.ParametersNotConformable);
            }

            return rightSide.LeftMultiply(leftSide);
        }
    }
}